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Symmetry- induced vortex-antivortex configurations in superconducting squares and triangles were 
predicted earlier; yet, they have not been resolved in experiment up to date. Namely, with vortex- 
antivortex states being highly unstable with respect to defects and temperature fluctuations, it is 
very unlikely that samples can be fabricated with the needed quality. Here we show how these 
drawbacks can be overcome by strategically placed nanoholes in the sample. As a result, (i) the 
actual shape of the sample becomes far less important, (ii) the stability of the vortex-antivortex 
configurations in general is substantially enhanced, and (iii) states comprising novel giant- antivortices 
(with higher winding numbers) become energetically favorable in perforated disks. In the analysis, 
we stress the potent of strong screening to destabilize the vortex-antivortex states. In turn, the 
screening-symmetry competition favors stabilization of new asymmetric ground states, which arise 
for small values of the effective Ginzburg-Landau parameter k* . 

PACS numbers: 74.20.De, 74.25.Dw, 74.25.Qt, 74.78.Na 



I. INTRODUCTION 

The question of vortex-antivortex (VAV) states in 
mesoscopic samples is one of the most intriguing ones in 
superconductivity studies of the last decade. Being sta- 
bilized in a homogeneous magnetic field, VAV states are 
counterintuitive and rather difficult to explain in usual 
terms. Nevertheless, by analyzing the solutions of the lin- 
earized Ginzburg-Landau (GL) theory, Chibotaru et al 1 
predicted their stability in fiat superconducting squares 
as a consequence of the symmetry of the sample. Namely, 
in linear theory, the geometry of the sample boundaries 
directly translates on the vortex states. As a conse- 
quence, for three fiuxons captured by the sample (i.e. 
vorticity L = 3), the vortex state with four vortices and 
one antivortex can be energetically more favorable than 
the triangular configuration of single vortices. Still, one 
should bear in mind that the linear approach of Refi is 
only valid at the superconducting/normal (S/N) bound- 
ary, where \£ is extremely small and the non-linear term 
in the GL theory becomes negligible. In a further step 3,4 , 
by taking into account the non-linearity of the first GL 
equation, the existence of VAV states was predicted to 
persist even away from the S/N phase boundary. The 
antivortex for L = 3 in a square of size 10£o x 10£o 
was shown to be stable in a temperature range of about 
AT = 0.3T C (with £o being the coherence length at 
T = 0, and T c the critical temperature). Up to now, 
such VAV configurations have not been realized in ex- 
periment. From a theoretical viewpoint, there are several 
reasons why these novel vortex states escaped experimen- 
tal observation: first, even a tiny defect at the boundary 
(< 1%) destroys the vortex-antivortex state 5 . Second, 
the vortices and the antivortex are all confined in a small 
area of typical size less than the coherence length £(T). 
As a consequence, this proximity results in a strong lo- 
cal suppression of the order parameter \I/ which makes 
the separate minima experimentally undistinguishable - 



the vortex configuration as a whole is very similar to a 
single giant-vortex. And third, imaging of the magnetic 
field profile is also shortcoming: the field generated by 
the sample is proportional to the supercurrent which is 
rather weak in the VAV region (being directly related to 

*)■ 

Therefore, for future experimental verification of the 
VAV states, it is essential: i) to increase the vortex- 
antivortex distance, allowing for measurements with re- 
alistic spatial resolution; ii) to obtain a larger contrast 
in the order parameter, and iii) the magnetic field, re- 
quired for e.g. Scanning Tunnelling and Hall-probe mi- 
croscopy measurements. However, to observe the 'pure' 
VAV nucleation, we should realize this in a homogeneous 
magnetic field, contrary to the suggestion of Ref. 6 where 
a magnetic dot was placed on top of the sample. It is 
already well-known that VAV states have low energy in 
an inhomogeneous magnetic field, such as the stray field 
of the magnetic dots 7 which may induce antivortices in 
its own. 

In our recent Letter—, we introduced the idea of strate- 
gically made perforations/holes/antidots in a square 
mesoscopic superconductor (see Fig. [1]) which act as 




FIG. 1: (Color online) Superconducting square with four 
nanoholes (i.e. antidots). 
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pinning centers for the vortices. The pinning force is 
supposed to pull the vortices further from the sample 
center and thus create a larger separation between the 
vortices and the central antivortex. In what follows, 
we will elaborate further on this concept, and demon- 
strate a number of its advantages. Contrary to Ref. 8 
where we considered only the square geometry, here we 
expand the study to different polygonal shapes of the 
sample and geometry of the pinning. In addition, while 
previous studiesi* 2 -^ only studied VAV configurations 
in extremely thin samples, we consider samples of fi- 
nite (though relatively small) thickness, by incorporating 
magnetic screening effects in the calculation. 

The paper is organized as follows. First the theoretical 
approach is formulated (Sec. II), comprising the full non- 
linear GL theory and the linearized GL theory (LGL). In 
Sec. Ill we present our results, and determine the op- 
timal pinning parameters to enhance the stability range 
of the vortex-antivortex configuration. The influence of 
imperfections and defects on the superconducting state 
is discussed in Sec. IV, the competition of different sym- 
metries in the sample is treated (Sec. V), e.g. a disk 
with a number of holes. Next we describe the influence 
of the non-linear term in GL theory and GL parameter 
ft, i.e. screening (Sec. VI) and some guidelines for a 
possible experimental observation (Sec. VII). Finally, we 
summarize our findings in Sec. VIII. 



II. THEORETICAL FORMALISM 

A. Ginzburg-Landau theory 

It is already well established that the Ginzburg-Landau 
(GL) formalism gives an accurate description of the su- 
perconducting state of low-T c superconductors, and is 
well suited to incorporate boundary effects in the treat- 
ment of the mesoscopic superconductors. The GL theory 
postulates the following expression for the free energy of 
the system (in dimensionless form): 



J 



dr 



-2|*| 2 + |*| 4 + 2|(-iV- A)¥|' 
+2k 2 (B - B appl f 



(1) 



Here all distances are measured in units of the coher- 
ence length f (T) = £o/\/l - T/T c . The vector potential 
A is expressed in units of cft/2e£, the magnetic field in 
units of H C 2 = ft/2e£ 2 , and the complex order parameter 
^ in units of = \/ — such that = 1 in the 
pure Meissner phase and = in the normal conduct- 
ing state (with a, (3 being the GL coefficients^). |\I/| 2 
represents the local Cooper-pair density, which we will 
abbreviate as CPD. T5 ap pi denotes the applied magnetic 
field, while B stands for the total magnetic field, ft is the 
ratio between penetration depth and coherence length, 



i.e. ft = A/£ and is assumed temperature independent in 
this model 

To find stable solutions, one has to find the wave func- 
tion and the vector potential which minimize the free 
energy functional. Two coupled non-linear differential 
equations, one for the order parameter and one for the 
vector potential, can be derived using variation analysis. 
We are interested in thin, but finite- thickness, samples in 
a perpendicular applied field, where we may neglect the 
variation of the magnetic field and order parameter over 
the thickness of the sample. Accordingly, we average the 
GL equations over the sample thicknes G 13 i 14 i 15 , and write 
them as 
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where j s denotes the supercurrent density, and ft* = 
ft 2 /(d/£) is the effective Ginzburg-Landau parameter, 
the consequence of averaging over sample thickness d. 
It's worth emphasizing that ft* is temperature dependent 
(e.g. through £(T) in the denominator). All considered 
samples in this paper, treated with the full non-linear GL 
theory, have the same thickness d = £o- Consequently ft* 
takes the temperature dependence ft*(T) = ft 2 /Vl — T. 
To address properly the influence of T, all results ob- 
tained by the full GL theory will be expressed in units of 
£, and H C 2 at zero temperature. 

While authors of Refsi* 2 -^ solved only the first GL 
equation, and restricted themselves to the limit of d — > 0, 
or extreme type II superconductors, we solve numerically 
both equations. The strength of the influence of the sec- 
ond GL equation is governed by ft* . When ft* ^> 1 (thus 
for d <C 1 and/or ft ^> 1), Eq. (j3]) can indeed be ne- 
glected. However, when ft* « 1, we found it can have a 
tremendous influence and provide us with some interest- 
ing and novel behavior of the superconducting state. 

After solving Eqs. ([2j[3j) , the vorticity L of the found 
stable vortex configuration is defined as the sum of the 
individual winding numbers of all vortices present in the 
sample. The winding number of a vortex is defined as the 
integer number of 2ir phase changes of the order param- 
eter when circling around the vortex core. Consequently, 
the vorticity of a vortex is defined as 1, and the vorticity 
of an antivortex as -1. The total vorticity of a partic- 
ular vortex configuration can also be determined at the 
sample boundary fl by the formula 
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The vortex state of small mesoscopic superconductors 
is strongly influenced by the imposed topological confine- 
ment. The shape of the sample boundaries is introduced 
in our calculation through the Neumann boundary con- 
dition, which sets the supercurrent perpendicular to the 
boundary equal to zero: 
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n- (-iV-A)| boundary = 0, (5) 

where n is the normal unit vector on the surface. In 
Refs. 1,9,10 a unitary transformation was used in combi- 
nation with a certain gauge for the vector potential such 
that for the linear GL equation A = on the sample 
boundary. This approach is restricted to simple geome- 
tries (i.e. polygons) and cannot be applied to our multi- 
ple connected samples that contain holes. 

Our numerical approach is based on a finite difference 
Gauss- Seidel relaxation for the time dependent GL equa- 
tions. The numerical technique we applied is the one 
introduced in Refjii with convenient iterative expansion 
of the non-linear term explained in Ref J^, to speed up 
convergence. We solve the GL equations on a uniform 
Cartesian grid typically with 128 x 128 points. 

B. Linearized approach 

Since the VAV state is a consequence of the symmetry 
of the sample boundary and because of the most effi- 
cient reflection of this symmetry into the solutions of the 
linearized Ginzburg-Landau (LGL) theory, the LGL for- 
malism is the ideal instrument to study the influence of 
the geometry on the VAV state. From a numerical point 
of view, the LGL approach is far less demanding, and 
can still provide us with the minimal requirements for 
the realization of the VAV states. 

The linear theory is exactly valid only on the S/N 
boundary, i.e. when the |\I/| 3 term can be neglected. 
Due to the weak superconducting order parameter, the 
magnetic field equals the applied one and the second GL 
equation can be discarded. At the same time, Eq. ([2]) 
becomes 

(iV + A) 2 ^ = o^, (6) 

where a now is equivalent to an eigenvalue. An important 
property of Eq. (j6j) is that its solution is independent of 
the size of the sample, as long as the applied flux is held 
constant, i.e. the LGL method is a scalable theory. 

To be able to compare the solutions of the LGL ob- 
tained for different geometries, one has to make sure that 
they have the same normalization. We use the following 
normalization of the wave function, y J v dV\^\ 2 = 1, 
or equivalently (|^| 2 ) = 1, where the integration is per- 
formed over the sample volume V. This matches the 
uniform solution of = 1 in the whole sample in the 
absence of a magnetic field. 

Our aim is to use the LGL model deeper inside the su- 
perconducting state, even though it is strictly not valid, 
simply as a limiting case of the full GL theory. For that 
purpose we still derive the generated magnetic field pro- 
file using Eq. (j3j). Although the LGL theory is not able 
to determine an absolute order of magnitude of the field 



(because ^ « 0), it gives the correct magnetic field up to 
a multiplication constant, and thus allows for comparison 
of field profiles for different sample geometries. 

We compute the vector potential A from the calculated 
distribution of supercurrents, for taken k* = 1. For this 
purpose, we solve two independent Poisson equations (for 
A x and A y ) using the Fast Fourier Transform. Magnetic 
field is then obtained as B = V x A. The relative mag- 
netic field profile found in this way turns out to be rather 
accurate, even when going deeper in the superconducting 
state (e.g. by decreasing temperature). When compared 
with a calculation based on the full non-linear GL theory, 
the main effect of the second equation is to enhance the 
magnetic field generated by the sample. 

As an advantage, to solve Eq. ([6]) is computationally 
not very demanding, and we use the numerical package 
COMSOL (formerly known as Femlab) for this purpose. 
This software package uses the finite element method to 
solve differential equations. Its profound precision (in 
solving linear equations) allows us to significantly in- 
crease the grid resolution, for instance, in the neighbor- 
hood of the VAV complex. The use of triangular finite 
elements proved as a more accurate treatment of the sam- 
ple's geometry than we could achieve with a rectangular 
grid, for an arbitrary sample geometry. 

Therefore, our linear approach is well suited for com- 
parison of the VAV state in different samples. We will 
use it to determine the optimal hole parameters and to 
study the influence of imperfections and defects. We will 
also treat the competition of different pinning symme- 
tries and samples with different geometries, like e.g. per- 
forated disks, in search for novel VAV states. 



III. ENGINEERING OF THE ARTIFICIAL 
PINNING SITES 

In what follows, we will study the optimal parameters 
of the artificial pinning centra, i.e. size and position of 
the holes introduced in the mesoscopic sample, using the 
linear GL theory. Later on, the influence of the non- 
linearity of the GL equations and the inclusion of the 
second GL equation will be systematically investigated. 

Our aim is to enhance the parameters that characterize 
the VAV state. These are: (i) the maximal Cooper-pair 
density (CPD) between the V and AV (\ip\max) because 
of the need of CPD imaging techniques (e.g. STM) for a 
sufficient contrast in the local density of states, (ii) the 
distance between V and AV (d vav ) which has to be larger 
than the spatial resolution of the measurements, and (iii) 
the magnetic field difference between V and AV (AB vav ). 
All the latter properties for a given VAV state are defined 
for the VAV pair with smallest separation d vav . 

The influence of the holes on the vortices is determined 
by their position and their size. How can we enhance 
e.g. the vortex- ant ivortex separation by manipulating 
these parameters? Holes are known to be pinning sites, 
i.e. they attract vortices. By placing four holes, we can 



4 




d h /D [%] 

FIG. 2: (Color online) Blue solid curve: The minimum size 
Wh of the hole for the vortices to be fully captured by the 
holes at distance dh- Above the blue line, vortices are con- 
tained in the holes. Below the blue line, vortices are outside 
and between the holes. Red dashed curve: The maximum 
CPD between the vortex and the antivortex when the holes 
are minimum sized. Green dotted curve: Dimensionless quan- 
tity AB vav , expressing the strength of the magnetic field dif- 
ference measured in the vortex and in the antivortex. Insets 
are contour plots of the phase of the wave function (blue/red 
corresponds to 0/2tt) for d h /D = 28% and w h /W = 18% 
(upper) and d h /D = 36% and w h /W = 14% (lower). The 
applied flux = 5.50o. 
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FIG. 3: (Color online) Blue solid curve: The minimum size 
Wh of the hole for the vortices to be fully captured by the 
holes at distance dh- Above the blue line, vortices are con- 
tained in the holes. Below the blue line, vortices are outside 
and between the holes. Red dashed curve: The maximum 
CPD between the vortex and the anti vortex when the holes 
are minimum sized. Green dotted curve: Dimensionless quan- 
tity AB V av, expressing the strength of the magnetic field dif- 
ference measured in the vortex and in the antivortex. Insets 
are contour plots of the CPD of the wave function (blue/red 
corresponds to low/high density - logarithmic color scale) for 
dh/D = 22% and w h /W = 24% (upper) and d h /D = 24% 
and u>h/W = 18% (lower). The applied flux <j) — 3.50o- 



pull the vortices, which surround the antivortex, away 
from the antivortex. The larger the hole, the stronger 
the attraction. The forces on the antivortex are exactly 
cancelled because of the symmetric position of the holes. 
However, the vortices also experience inward forces, i.e. 
they are attracted by the antivortex and additionally the 
Meissner current also compresses the vortices to the in- 
side of the sample. The vortices eventually will find an 
equilibrium position, in which the inward and outward 
forces cancel exactly. This position depends of the dis- 
tance of the hole to the center and of the size of this 
hole. 

To understand the influence of the hole parameters 
(size and position) we introduce a criterium which we can 
use as a guide to determine the optimal hole parameters. 
This criterium is: What is the minimum hole size for a 
given hole position, such that the vortex is still captured 
by the hole? It seems that when holes are made larger, 
they attract the vortices more strongly and eventually 
the vortex will be captured by the hole, when the holes 
are placed not too close to the boundary of the sample. 
We show our results in Fig. [2] and Fig. [3] for the square 
and the triangle geometry, respectively. The insets show 
the setup and the definition of the hole parameters. The 
characteristic variables dh(= d vav ), \ip\max an d AB vav 
are also shown in Figs. [2]and[3j 

We note that the curves of 1^1 max and AB vav are iso- 
morphous. We expect this, because max ex \i(j\ 2 ex 



j s ex AB cx AB vav . The maximum 1^1 max anc ^ as a 
consequence, AB vav in a triangle appear to be 40 % less 
than for the square. 

From these figures, we can conclude that, for both the 
square and the triangle shaped samples, moving the holes 
farther from the center implies enlarging them, in order 
to keep the vortices trapped by the holes. From a cer- 
tain distance dh, the needed hole size Wh diverges, i.e. 
from a certain distance of the holes, they can never be 
made large enough to trap the vortices that surround the 
antivortex. For both the square and the triangle, this dis- 
tance dh/D is about 30%. While for the square increas- 
ing the distance always implies an increase of max and 
AB vav , this is not true for the triangle: from a distance 
dh/D « 23% these observables start degrading. 

However these results apply to a specific shape of the 
holes, and in principle it is not fair to compare the perfo- 
rated square with the perforated triangle, since the holes 
have different shape. For the square, the trapped vortex 
is located in the corner of the square hole; for the triangle 
this vortex is located in the middle of a side. Thus, as 
an alternative approach we reinvestigated both samples, 
now with circular holes, and noticed that qualitatively 
these results do not differ from the ones with the other 
shaped holes. Consequently, we conclude that the gen- 
eral conclusions for square and triangle are independent 
of the shape of the holes. 

When, in a realistic situation, one is interested in cap- 
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FIG. 4: (Color online) VAV distance (d vav ) and the maximal 
value of the CPD between the VAV (\ip\max) versus the hole 
size Wh for a square with square holes at position dh/D = 
28%. Applied flux = 50 o . 



turing the vortices by the holes, one has to chose a some- 
what bigger hole size than the minimal one since small 
defects will distort the VAV state and drive a VAV pair 
closer to each other, hereby the vortex leaving the hole. 

For holes at a distance of dh/D = 28%, the influence 
of the hole size is depicted in Fig. [U However a larger 
hole size will of course not increase the VAV distance, it 
will increase \^\^ ax which is important for the needed 
experimental resolution. 



IV. INFLUENCE OF 
IMPERFECTIONS /DEFECTS 

By introducing nano-engineered holes, we enhanced 
the observability of the VAV state. But how about its 
stability against defects? 

It is generally known that the VAV state in a plain 
square is very vulnerable to defects positioned at the 
edge of the sample 5 . We will show how this situation 
dramatically changes when fourfold symmetric holes are 
introduced into the sample. We will consider defects as 
slight modifications of the geometry using the linear GL 
theory. We observe that when the linear theory predicts 
that the vortex- antivortex pair annihilates, it will surely 
be so in the non-linear theory, which justifies the use of 
the LGL theory, in order to determine the minimal re- 
quirements for the stability of the VAV state. 



A. A defect at the edge 

As a starting point for comparing the perforated sam- 
ple with the plain sample we will study the influence of 
defects at the edge of the sample. It is generally known 
that a small defect at the boundary will destroy the VAV 
state in a plain square 5 -. In the present study we restrict 




p d =o 




□ 




□ 


□ 




10 



15 



w d /W [%] 



FIG. 5: (Color online) The vortex-antivortex distance vs de- 
fect size for a plain (dark grey) and a perforated (light grey) 
square, with a square defect at the middle of one edge, i.e. 
Pd = 0. The dashed curve represents an indentation, the 
solid curve a bulge defect. The applied flux for the plain 
square <j> — 5.50o and for the perforated square <f> — 4.O0o 



ourselves to defects that are indentations and bulges at 
the surface of the sample. We will present here the re- 
sults for a square for two different defects positioned at 
the edge of the sample. The defect under study is placed 
at the center of the edge of the square and is taken to be 
a square itself with side Wd- We considered a small bulge 
and a small indentation. It is known that such defects 
may influence the penetration and expulsion of vortices 
as was demonstrated experimentally by A. K. Geim et. 
al in Ref. 17 . For theoretical studies we refer to Refs. 18 
and 19 . 

We focus on the VAV separation d vav to compare the 
influence of defect size between the perforated square and 
the plain square. When this distance equals zero, the 
VAV pair annihilates and the VAV configuration disap- 
pears. The results for the plain square are shown by the 
two curves at the lower left region of Fig. [5j The results 
for the perforated square are shown by the two curves 
at the upper right region. The hole size and position 
used for the perforated square are Wh/W = 20% and 
dh/D = 28%. For the plain square we applied a flux 
(j) = 5.50o, for the perforated square <\> = 4.O0o- We had 
to apply different fluxes since the stability range of the 
L=3 state in both samples do not overlap. (The intro- 
duction of holes shifts the phase diagrams to lower flux.) 

For a plain square with a bulge defect, we conclude that 
the VAV configuration exists until a size Wd/W ~ 5%. 
The VAV state is even more sensitive to an indentation: 
a size of ~ 1% is sufficient for the disappearance of the 
VAV. However, we found that the perforated square is 
much more resistant to defects. Bulge and indentation 
defects, up to 14% (comparable to the size of a nano- 
hole), now coexist with a VAV configuration. Bulge and 
indentation defects have a similar effect on the VAV sta- 
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FIG. 6: (Color online) The influence of a defect on the vav 
distance (d vav ), CPD between vortex and antivortex ^ax) 
and the magnetic field difference between vortex and antivor- 
tex (AB V av)- The hole parameters are Wh/W = 20% and 
dh/D = 28% and the applied flux = 40o-- 
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FIG. 7: (Color online) Contourplots of the phase of the order 
parameter in a perforated square with several kinds of defects. 
The hole parameters are: w h /W = 20%, d h /D = 22% and 
applied flux <j> — 4.10o- The defects include: a) a diagonally 
displaced hole, b) a horizontally displaced hole, c) a bigger 
hole and d) a different shaped hole. 



bility, in contrast to the plain square case. We conclude 
that the 4-fold symmetrically placed pinning centers are 
much more efficient in imposing their symmetry than the 
outer boundary of the sample. As a consequence, small 
defects at the edge have little effect and only distort the 
VAV molecule 

One may argue that this effect is possibly only since the 
defect is located exactly in the center of the edge, there- 
fore imposing mirror symmetry. Because of this symme- 
try, the antivortex is prohibited to chose one of the two 
vortices of the holes to annihilate with and it stays on 
the mirror line. 

To turn off this stabilization effect due to symmetry, 
we also investigated defects which are displaced from the 
middle of the edge. We studied the influence of the posi- 
tion and size of an edge defect on the characteristic pa- 
rameters d vav , l^lmaa? an d &B vav . The result is shown 
in Fig. [6] for a perforated square with hole parameters 
Wh/W = 20% and dh/D = 28% under an applied flux of 
4(/>o- On the left side, the influence of the defect size is 
depicted for both a bulge and an indentation defect, po- 
sitioned at a distance Pd/W = 3% from the edge center. 
On the right side, the influence of the position of a defect 
of size Wd/W = 7% is shown. 

Concerning the influence of the size of the defect, we 
notice that bulge and indentation defect act similarly: In- 
crease of the defect size implies degradation of all three 
observables l^lma^ d vav and AB vav . However, for the 
position of the defect, we predict different behavior for 
bulge and indentation. For a bulge defect, we notice that 
the stability of the VAV configuration decreases as the 
defect moves further from the center of a side, while for 
the indentation we observe a decrease followed by an in- 
crease. However, for both types of defects, the VAV con- 



figuration survives best when the defect is centered at 
an edge. We also conclude that a bulge defect is less 
destructive than an indentation, unless it is placed near 
the center of an edge. The antivortex acts like being 
attracted to the indentation and repelled by the bulge 
defect. 



B. Other kinds of defects 

The VAV stability against edge defects is strongly im- 
proved by the introduction of the fourfold symmetrically 
placed holes. However, in actual experiments these holes 
themselves can contain defects or imperfections such as 
non-uniform sized holes, holes which are a little displaced 
with respect to their high symmetry position, or holes 
with a slightly different shape. 

In Fig. [71 four examples are shown. In all pictures 
the 'normal' holes can be described by the parameters 
Wh/W = 20%, dh/D = 22% and we apply a magnetic 
flux <p = 4.10o- The plots in Figs. [71 (a) and (b) il- 
lustrate the effect of imperfect positioned holes. In Fig. 
(a), a diagonal displacement of the upper right hole 
over a distance d/W = 6% is shown. The VAV survives 
such displacements in the range of -4% to +12%. For a 
horizontal displacement, as illustrated in b), the range is 
smaller: from -3.5% to +3.5%. The fact that the VAV 
configuration is more stable for a diagonal displacement, 
we attribute to the existence of mirror symmetry along 
one diagonal, forcing the AV on this diagonal. 

In FigEl (c) the upper right hole has a somewhat larger 
size than the other holes. Its sides are 5% larger. This 
kind of hole size defect does not destroy the VAV state, 



FIG. 8: (Color online) Contourplot of the phase of the order 
parameter for a square sample with three circular holes. The 
state with vorticity 2 is shown. The symmetric positioning of 
the holes induces an antivortex in the center. 



when in the range of -8% to +9%. 

In Fig. (d) the upper right hole is circular. It has 
the same area as the square holes and is centered like the 
square holes. The VAV survives, which illustrates that 
it's not the exact shape of the holes which matters, but 
rather its area. This implies that non perfect holes, i.e. 
holes with defects, will not destroy the VAV state as long 
as the imperfection is not too large. 



V. COMPETING SYMMETRIES 

A. Superconducting samples with polygonal 
pinning 

Several symmetries compete with each other when the 
VAV state nucleates. First we have the vortex- vortex 
interaction through which the vortices will try to form 
the Abrikosov triangle lattice. Second there is the sample 
boundary which will impose its own symmetry, and third 
there are the pinning centers (i.e. the holes) which will 
also try to impose their own symmetry. 

For small samples of the order of several coherence 
lengths, like the ones we studied up to now, the symmetry 
of the sample boundary opposes the Abrokisov lattice. In 
this subsection we'll point out that the symmetry of the 
pinning sites is even more strongly dominating than the 
symmetry of the boundary. 

For example we show the VAV state with vorticity 
L = 3 — 1 = 2 in a square with three holes in Fig. [8j 
Although the outer edge has a square geometry it is the 
triangular arrangement of the circular holes which im- 
poses its symmetry on the superconducting wave func- 
tion in the center of the square. This leads to 3 vortices 
trapped in the holes and a single antivortex in the middle 
of the sample. 



FIG. 9: The parameters R, p, r, d and used to characterize 
the geometry of the disk with N holes. 




FIG. 10: (Color online) Disk with 20 holes, in the state with 
vorticity L=19, consisting of 20 vortices and 1 antivortex in 
the center. Left: Logarithm of CPD Contour plot (Red (Dark 
grey) /Yellow (Light grey) = High/Low density) Right: Phase 
Contour plot (Red (Dark grey)/Yellow (Light grey) = 0/27r) 



B. Superconducting disk with many holes 

The arrangement of holes, instead of the sample outer 
boundary, seems to be the geometry element which is 
much more effective in imposing its geometry on the wave 
function. For this reason, we can as well use circular sym- 
metric disks, perforated by N symmetrically placed holes, 
and expect an N-fold symmetric CPD. A large variety of 
symmetry induced antivortex configurations can be cre- 
ated this way. Following this approach we found that it 
was possible to even create giant anti-vortices with vor- 
ticity L up to —7. 

The parameters which define the geometry are de- 
picted in Fig. [9j d is the distance between holes, the 
angle between holes which is set equal to 2tt/N, as to im- 
pose the N-fold symmetry. There is also p, the distance 
of the holes' center to the disks center. R is the radius of 
the disk and r is the radius of the holes. Imposing that 
the distance between the holes equals the hole's diameter 
is equivalent to the condition r = (l/2)psin(#/2). 

To illustrate this conjecture, we give some examples. 
In Fig. [10] the CPD and the phase plot of a disk with 
20 holes is shown. The parameters are r/R = 4.3%, 
p/R = 55%, so that the distance between the holes 
d = 2r. The applied flux is <p = 26cV The vorticity 




FIG. 11: (Color online) Disk with 20 holes, showing 4 antivor- 
tices in the center. Left: Logarithm of CPD Contour plot 
(Red (Dark grey) /Yellow (Light grey) = High/Low density) 
Right: Phase Contour plot (Red (Dark grey) /Yellow (Light 
grey)=0/27r) 



L = 19, is one less the symmetry of the holes and con- 
sequently an additional vortex-antivortex pair is created. 
The alternative solution to obey the symmetry would be 
to create a L = 19 giant vortex, but it turns out that the 
creation of one VAV pair is energetically more favorable. 

There's clearly a competition between the formation 
of a giant vortex - the merging of several vortices in one 
point - and the formation of a giant antivortex - consist- 
ing of the creation of vortex-antivortex pairs and then 
merging of all anti- vortices in a singular point. 

For the same disk, with an applied flux of <j) = 210o and 
slightly modified parameters p/R = 0.5 and r/R = 4.73% 
the CPD and the corresponding phase plot is shown in 
Fig. [TTJ The phase plot indicates a giant antivortex in 
the center with winding number -4. We found that the 
giant antivortex is highly unstable with respect to de- 
fects. Therefore, in experiment, it will be very difficult, 
to observe these giant antivortices, and most likely they 
will fall apart in 4 single antivortices. The reasons for 
this are: i) a sample with a perfect symmetry and in the 
absence of defects will be difficult to manufacture. In a 
first stadium, the imperfections will cause the giant AV 
to disassociate into separate antivortices and for larger 
imperfections it will cause the VAV pairs to annihilate, 
ii) The extreme low CPD in the center due to the densely 
packed (anti-)vortices which is almost impossible to dis- 
tinguish between the separate minima. 

To generalize the concept of giant-antivortices and 
their appearance in disks with N symmetrically placed 
holes, we investigated the relation between the symme- 
try order N and the vorticity L. The result is the L-N 
phase diagram shown in Fig. [12] where we used the fol- 
lowing parameters: p/R = 0.55, r = 2d and <p = 1.2N<pQ. 

The anti- vorticity of the state is indicated by the differ- 
ent colors. When L < N and there are no anti- vortices, 
a giant vortex of vorticity L is located at the disk center, 
since this is the only solution which is able to obey the 
symmetry. However, because of the finite grid resolution, 
these vortices split up in separate ones, analogous to the 
fate of the giant antivortex described before. 

For vorticities L > nN (where n is the largest inte- 
ger obeying the condition) another property is observed: 




5 10 
Vorticity (L) 

FIG. 12: (Color online) For each vorticity and for each num- 
ber of holes N, the number of anti- vortices is shown. Due to 
the symmetry the antivortices combine to a giant antivortex. 
Left of the thick line no antivortices are observed. 



every hole pins n vortices, which are not necessarily cap- 
tured inside the holes, but clearly belong (i.e. are at- 
tached) to the specific hole. 

The L-N phase diagram is not uniquely determined by 
N. The choice between giant- or anti-vortices is strongly 
dependent on the exact choice of the geometry (i.e. of 
r, p, d) which strongly affect the free energy, and last 
but not least of the applied flux, because of the Meissner 
current compressing the (anti-) vortices inward. 



VI. INFLUENCE OF THE MAGNETIC 
SCREENING 

The non-linearity of the GL equations and the coupling 
to the magnetic field will now be taken into account. This 
means that from now on, we will use the full non-linear 
GL theory. We will focus only on the square geometry, 
but all conclusions can be extrapolated to other geome- 
tries/symmetries. We remind you that the thickness of 
the sample throughout this paper is chosen to be l£o- 

The VAV state is now reached by, for instance, in- 
creasing the temperature. In this case the multivortex 
state will, through a second order transition, transform 
into the symmetric VAV state. In between these two 
phases, a new state arises: the asymmetric VAV state. 
This is a stable ground state configuration, consisting 
of several vortices and one antivortex, but they are not 
positioned symmetrically. The area in a <j) — T phase dia- 
gram, where these asymmetric symmetry-induced states 
are the ground state configuration, increases when k de- 
creases, and their existence is therefore a consequence of 
the second GL equation. 

As a definition for the lower temperature of the asym- 
metric VAV region, we take the start of the nucleation of 
the vortex-antivortex pair. See inset B in Fig. [16] for an 
illustration of the CPD at this point. In our simulations 
we noticed a sharp and sudden drop of the CPD in the 



9 



0.95 
0.90 
0.85 
0.801 
|_° 0.75 
0.70 
0.65 
0.60 
0.55 



Normal State 




mm 



I I Multi vortex 
W\ Asymm. VAV 
□ Symm. VAV 



o 



5 10 

W/W [%] 



15 



FIG. 13: (Color online) Transition temperatures T vav , 
Tvav-symm, ^s/n vs. the size of the hole for a square of sizes 
10£o x 10£o and an applied field of <j)/(j)Q — 5.0. The corners 
of the holes are at a distance dh/D = 25%. On the left side 
of the white line, the vortices are captured by the holes, on 
the right side they are outside the holes, k = 1.5 



point where later the vortex and antivortex are separat- 
ing. At this point where the phase change around this 
pixel is still zero, we define the beginning of the VAV nu- 
cleation when the CPD drops with a factor < 2OOO0o/^c- 
We say that a VAV state is symmetric, when a, our mea- 
sure of asymmetry is smaller than 0.05, which is defined 



as a 



Z)min \xi 

i 3& 



i/few). 



Here, R is the ro- 



tation operator, which rotates a vector over 90 degrees, 
x$ are the positions of (anti-)vortices (when they're con- 
tained in holes, the position is taken to be in the inner 
corner pixel of the hole). The origin of the axes is chosen 
in the middle of all vortices, i.e. x i = 0- 

i 

In the non-linear theory we will now review the role 
played by the hole size. In Fig. [13] a phase diagram shows 
the dependence of the VAV stability temperature interval 
of the hole size. In this figure we used n = 1.5 and the to- 
tal applied flux <\> equals 5.O0o- The position of the inner 
corner of the holes was fixed at dh/D = 25%. Although 
increasing the hole size, enhances the characteristic ob- 
servables |^|^ a;E and AB vav and as well stabilizes against 
all kinds of defects, we notice one backdraw of large holes: 
the VAV temperature interval shrinks. One can also see 
that the VAV temperature interval reaches a maximum 
size at about Wh/W = 5% while the temperature range of 
asymmetric VAV states remains almost unaltered before 
this maximum. 

(j) — T phase diagrams are shown in Fig. (TUfor a plain 
square and in Fig. [15] for a perforated square, both of 
size 10£o x 10£o- Both figures contain two superimposed 
phase diagrams, one for n = oo (indicated by grey) and 
one for k = 1 (shown in red (dark grey)). The diagrams 
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FIG. 14: (Color online) cj>— T phase diagram for a plain square, 
illustrating the ground state regime of the L=3 state. Grey 
represents k = oo, red (dark grey) k = 1. 





1 1 1 1 


i i 1 i i i i 1 i i 

- 1=4 


1 










\ 

\ 

\ 

\ 


^^^^^^ 




\ 

\ 

\ 

\ 


1 1 Multivortex 






7771 Asymm. VAV - 




i 1 i i 


■ Symm. VAV ■ 

\ i \i i i i i i i i 



0.9 
0.8 
0.7 
3 0.6 
0.5 
0.4 
0.3 



FIG. 15: (Color online) <j> — T phase diagram for a square 
containing four symmetrically placed nanoholes, illustrating 
the ground state regime of the L=3 state. Grey represents 
k — oo, red (dark grey) k = 1. Hole parameters: Wh/W — 
12.5% and d h /D = 25%. 



illustrate the ground state region of the L=3 state, so 
only the neighboring L=2 and L=4 states are taken into 
account. 

Comparing these two phase diagrams with and without 
holes, one can clearly see that the introduction of holes of 
size Wh/W = 12.5% at the position dh/D = 25% causes 
the total VAV region to shrink for a sample with n = oo 
but to expand for a sample with k = 1. Nonetheless, 
the region of asymmetric VAV states shrinks for both k. 
Note however that this is not true for all sizes. See e.g. 
Fig. [13] where the asymmetric region initially undergoes 
a subtle increase. 

For both squares a decrease of k leads to: i) a broaden- 
ing of the temperature regime with an asymmetric vortex 
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FIG. 16: (Color online) k-T phase diagram showing the 
regions where symmetric (cyan/dark grey), asymmetric 
(shaded) and no antivortex (yellow/light grey) states can be 
found. Plain square 10£o x 10£o- Flux — 6. Insets 

show contour plots of the logarithm of the CPD for the three 
regions. 



state, and ii) it causes a shift of the vorticity to higher 
fields which is large for low temperatures and smoothly 
decreases to zero, at the S/N boundary. The reason 
for the latter is that the larger currents (j s cx |^| 2 ex 
1 — T/T c ) generate a stronger magnetic field and thus are 
more effective in expelling and concentrating the mag- 
netic flux. This shift for decreasing k can be understood 
in the limit to type I superconductivity, where L=0 (i.e. 
the Meissner state) is the only stable state. The S/N 
boundary of both n coincide, since at the S/N boundary 
the second equation does not have an influence. 

Fig. [16] depicts the ^-dependence of the temperature- 
interval wherein the (asymmetric) VAV state is stable 
in a plain square. It seems that the T-interval of the 
asymmetric VAV state grows for decreasing k. Also it 
is shown that the decrease of n disfavours the vortex- 
antivortex nucleation, which is opposite to the findings 
of Ref. 12 where this conclusion was made for a mesoscopic 
type-I triangle. 

For low n* , the supercurrents generate a magnetic field 
which contributes to the total magnetic field and destroys 
its homogeneity. Obviously, the shape of the magnetic 
field shows strong similarities with that of the Cooper 
pair density, which could be the reason for the stabiliza- 
tion of asymmetric VAV states for low k*. 

From this point of view we can also interpret the broad- 
ening of the stability region of asymmetric VAV states for 
small n. The decrease of n acts similarly as having a de- 
fect, but one without preferential spatial direction. It 
only slows the nucleation/ annihilation of the VAV pair. 

Since the vortex- antivortex state already is a state 
which is rather unstable and sensitive to all kinds of de- 
fects, it is normal that this subtle equilibrium of the co- 
existence of vortex and antivortex disappears and that 
the VAV pair annihilates. 




FIG. 17: (Color online) Contourplot of the quantity 
ABvav = (By — Bav)/ B app i as function of k and T for 
cf) = 5.8, T = 0.84 and k = 0.7 for a 10£ x 10£ supercon- 
ducting square with holes. Hole parameters: Wh/W = 12.5% 
and d h /D = 25%. 



From a theoretical viewpoint, there are several ways to 
make the second order transition from the highly sym- 
metric VAV state to the multivortex state: by decreas- 
ing temperature, the magnetic field or n. In all these 
three scenarios the same transition takes place qualita- 
tively: The antivortex moves towards one of the vortices, 
they approach and eventually annihilate. This clearly is 
a manifestation of spontaneous symmetry breaking. 



VII. HOW TO DETECT THE VAV STATE? 

Essentially, we see two ways to proof experimentally 
the existence of stable antivortex states: through mag- 
netic field imaging (e.g. using a Hall-probe) and through 
Cooper pair density imaging (e.g. using Scanning Tun- 
neling Microscopy). In this section we will highlight the 
advantages and disadventages of both approaches. 

An antivortex is characterized by the following prop- 
erties: i) the CPD in the center is suppressed to be ex- 
actly zero, and ii) its supercurrents circulate opposite to 
the one of the vortices. Unfortunately, the first property 
also applies to a conventional vortex, so that it cannot 
be used to discriminate between vortex and antivortex. 
This means that the magnetic field, generated by the 
supercurrents, is the only observable parameter able to 
distinguish between vortex and antivortex. 

In our study of superconductors in a homogeneous 
field, antivortices always seemed to appear surrounded by 
vortices. However, (multi)vortex currents always gener- 
ate in the center of the sample an anti- vortexlike current, 
in a superconducting cog wheels motion of the conden- 
sate. This current masks the current profile of a possible 
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FIG. 18: (Color online) Contourplot of the quantity A5 = 
(B - B app i)/B app i for = 5.8, T = 0.84 and k = 0.7 in a 
10£o x 10£o superconducting square with holes in the L=4-l 
state. Hole parameters: Wh/W = 12.5% and dh/D = 25%. 




FIG. 19: (Color online) The distribution of the supercurrent 
in a perforated square. The black arrows indicate the di- 
rection of the current, the background greyscale map depicts 
the magnitude of the local supercurrent density (white/black 
represents low/high). The main characteristics of the current 
flow pattern are indicated by the blue and red circles. Left: 
L = 4 state with applied flux cj) — 70o- Right: L = 4 — 1 = 3 
(VAV) state with applied flux cj) — 5.50o- 




FIG. 20: (Color online) Magnetic field profiles (diagonal and 
vertical cut) for the L=4-l=3 and L=4 states in a perforated 
square. AB = B — B app i- 



netic field profile. At low temperature the magnetic field 
profile should clearly exhibit the penetration of the field 
through the three vortices residing in three of the holes, 
while at a higher temperature the magnetic field should 
be penetrating the four holes equally, since a VAV pair 
is then created, one vortex occupying the empty hole. 

For a more direct observation of the antivortex, we 
suggest taking a sample with low effective ft* where a 
gradual second order nucleation of the vortex- antivortex 
pair is found. Then, using our findings from Sec. VI, 
subsequent CPD images for slowly increasing/decreasing 
temperature or magnetic field should demonstrate clear 
evidence for VAV-nucleation/annihilation. In this sce- 
nario, the hole parameters have to be chosen to max- 
imize the stability region of asymmetric VAV states in 
the vortex-phase diagram, following the guidelines from 
the preceding section. 



VIII. CONCLUSIONS 



antivortex which could be present since it has the same 
direction. As a consequence, we can not get evidence for 
the existence of an antivortex by looking in a qualitative 
way at the magnetic field profile of the sample. As an 
example, there is no qualitative difference of the super- 
current and the magnetic field profile between the L=4 
and L=3 VAV state, as is also shown in Figs. [T9landl20l 
Of course there is a quantitative difference, which could 
be exploited. 

Nonetheless, one can make use of qualitative differ- 
ences of both the CPD and the magnetic field profile, 
at different temperatures. Here the challenge is to con- 
trol the vorticity of the sample and to assure that it stays 
constant while measuring and sweeping the temperature. 
This way one can either measure the CPD or the mag- 



The existence of geometry induced antivortices in the 
presence of a homogeneous magnetic field has been pre- 
dicted theoretically several years age 1 . Up to now, there 
does not exist any experiment with the observation of it. 
With numerical simulations, using both the linear and 
full non-linear Ginzburg-Landau theory, we investigated 
how to engineer superconducting samples to stabilize and 
enhance the antivortex. 

We propose how to engineer the superconducting sam- 
ple but without taking away the conceptual novelty of the 
nucleation of the vortex- antivortex pair in a homogeneous 
magnetic field, as opposed to the idea of placing e.g. a 
magnetic dot on top of the sample. We pursued the idea 
to introduce holes which will act as pinning centers, and 
in doing so, pull the vortices away from the ant i- vortex 
and additionally to provide a strong immunity against 
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imperfections and defects. First we elaborated on the 
size and position of these holes. We determined optimal 
parameters for the square and triangle geometry. For in- 
stance, in a square geometry, we managed to enlarge the 
separation between vortex and antivortex with a factor 
4, compared to the case without holes. 

We investigated the influence of several kinds of geo- 
metrical defects on the VAV state. For all the imper- 
fections we studied, we found that the holes cause a sub- 
stantial increase of the stability of the VAV configuration 
with respect to defects. With the technology which is 
nowadays available it is possible to manufacture samples 
with the desired precision to proof the existence of VAV 
states. 

The geometry induced antivortices are known to be a 
consequence of the symmetry of the sample. Therefore, 
we focussed on the competition of different sources of 
symmetry in the mesoscopic superconductor. We con- 
cluded that the pinning centers are by far the most ef- 
ficient in imposing their symmetry, while the role of the 
outer boundary of the sample has a less determining 
role. Thinking further on this line, we studied circular 
disks, perforated by a number of symmetrically placed 
holes, and found giant antivortices up to a vorticity of - 
7. However, we found that this spectacular configuration, 
is highly sensitive to imperfections. 

The effect of the non-linearity of the first GL equations 
and the magnetic screening represented by the second GL 
equation was critically examined, allowing us to investi- 
gate the influence of temperature and non-zero thickness 
of the sample on the VAV state. 



We constructed phase diagrams for different values of 
k for the perforated and the plain square system show- 
ing the stability region of the VAV state in the <\> — T 
plane. For low value of k the introduction of holes en- 
larges the stability region. We also found the remarkable 
appearance of asymmetric VAV states, stable in a wide 
region of the phase diagram. Such asymmetric states are 
counterintuitive since the VAV state is known to be a 
consequence of symmetry and yet manifests in an asym- 
metric way because of the important non-linearity of the 
GL equations. Nevertheless, the existence of asymmetric 
VAV states can be the key property for an experiment 
to proof the existence of VAV states. We revisited our 
discussion about the hole size in this new context and 
found out that the size of the holes has a large influence 
on the stability region of the asymmetric states. 

We showed that a small value of n* disfavors the 
vortex-antivortex state in all investigated geometries 
(square, perforated square, perforated triangle), indicat- 
ing that the findings of Ref. 12 where a vortex-antivortex 
state in a type-I triangle is predicted with large VAV- 
separation cannot be correct. We mainly concentrated 
our discussion on square samples but we believe that the 
main conclusions are also valid for triangle samples. 
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